Kinetically driven helix formation during the homopolymer collapse process 
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Using Langevin simulations, we find that simple 'generic' bead-and-spring homopolymer chains in 
a sufficiently bad solvent spontaneously develop helical order during the process of collapsing from 
an initially stretched conformation. The helix formation is initiated by the unstable modes of the 
straight chain, which drive the system towards a long-lived metastable transient state. The effect is 
most pronounced if hydrodynamic interactions are screened. 

PACS numbers: 36.20.Ey,82.37.-j,87.15.He 



Helices are among the most common structures in na- 
ture. They are simple shapes which can be stabilized 
rather easily with suitable potentials [J, [3, [3] . a-helices 
in proteins have the favorable property of being highly 
"designable" in the sense that they are robust to mu- 
tations and thermodynamically stable Hence he- 
lices are interesting structures from a static point of 
view. In this Rapid Communication, we report on a 
dynamic nonequilibrium phenomenon where they also 
appear spontaneously, for purely kinetic reasons, with- 
out special potentials or local constraints In com- 
puter simulations, we find that simple, "generic" polymers 
in suitable environments recoil into well-ordered helices 
when released from an initially stretched configuration. 
This happens despite the fact that the ground state of the 
polymers is definitely non-helical. The helices emerge as 
metastable transient states with a finite, temperature- 
dependent lifetime, which may become quite large. 

The dynamics of polymer collapse has been studied in 
generic model chains by a number of authors 0, 0, 
Our work differs from these studies in two aspects. First, 
we study quenches to very low temperatures, about 1/100 
of the Theta temperature. Second, we do not start from 
an equilibrated high-temperature conformation (a coil), 
but from a stretched chain. Hence we study "mechanical" 
quenches rather than temperature quenches. Such rapid 
mechanical quenches are realized, e.g., if the chain is at- 
tached at both ends to objects that move apart, by links 
that are weak enough that they eventually break up. 

We have observed helix formation for various differ- 
ent choices of potentials. In this Rapid Communication, 
we restrict ourselves to the simple case of Lennard-Jones 
beads connected by harmonic springs for clarity of pre- 
sentation. Beads that are not direct neighbors on the 
chain interact via a truncated Lennard-Jones potential 
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connected by bonds subject to the spring potential 

V bond (r) = a (r - r a f (2) 

with a = 100e/er 2 and r = 0.85er. The theta tempera- 
ture for these chains is of the order ksTg ~ e [hJ. We 
investigate the kinetics of chain collapse at temperatures 
that are 10-100 times lower, starting from an initially 
straight chain. As dynamical model for the motion of 
the beads we use Langevin dynamics with and without 
hydrodynamics. Inertia effects are taken to be negligi- 
ble [lH. The natural units of our simulation are defined 
in terms of the bead size a, the Lennard-Jones energy e, 
and the friction coefficient ( (see below). Based on these 
quantities, the time unit is r = (<J 2 /e. 

We begin with discussing the collapse dynamics with- 
out hydrodynamics. The solvent surrounding the chain is 
then effectively replaced by the friction £ and a Gaussian 
distributed stochastic force m acting on the monomer i, 
which fulfills the conditions [12J < fji >= and 

< mA^VjAt') >= 2 C k B T 8. l3 8 aP 5{t - t'), (3) 

with monomer indices i, j = 1 . . . N, cartesian directions 
a, j3 G {x, y, z}, and t, t' two given times. The equations 
of motion are 



Cn = fi + ffi. 



(4) 



where the constant Co is chosen such that the potential is 
continuous everywhere. Adjacent beads in the chain are 



They are integrated using an Euler algorithm with the 
time step A t = 5 x 10~ 4 r. The stochastic noise fji was 
implemented by picking random numbers with uniform 
distribution, but the correct mean and variance ([3]) at 
every time step [13]. 

Figure. Q] shows examples of a collapsing chain at 
moderate and very low temperature (T = 0.5e/fcs and 
T = 0.01e//cB, respectively). At moderate temperature, 
a variant of the well-known pearl-necklace scenario 0, 0] 
is recovered: The collapse starts with the formation of 
small globules at both ends of the chain, which subse- 
quently grow in size and finally merge into a single glob- 
ule. At low temperature, the scenario is entirely different: 
Instead of disordered globules, well-ordered helices ap- 
pear at the ends of the chain. This process is symmetry- 
breaking, since the model potentials do not impose any 
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FIG. 1: Configurations of an initially stretched chain with 
length TV = 200 monomers at different times t for the tem- 
perature T — 0.5e/fcs (top) and T = O.Ole/fcs (bottom). 



chirality. The helical structures then propagate along 
the chain, until they meet in the middle. The final chain 
structure is an almost perfect helix with possibly one de- 
fect in the middle, depending on whether the chiralities 
of the two merging helices match each other. It is char- 
acterized by the radius R ~ 0.5er, the pitch D ~ la, and 
has roughly n ~ 4 monomers per turn. 

Simple energy considerations show that these helices 
are by no means the ground state of the chain. The 
lowest energy state of an infinitely long chain, where sur- 
face effects can be neglected, is most likely a structure 
with close-packed parallel strands. A straightforward en- 
ergy minimization for this structure yields the energy per 
monomer U tot /N = — 7.75e. In contrast, the lowest pos- 
sible energy per monomer of an infinitely long helix is 
Utot/N = — 3.26e, which is clearly higher. This is of 
course not surprising, given that helices have a large, en- 
ergetically unfavorable polymer-solvent interface. They 
cannot possibly be ground states at infinite chain length. 
The question remains whether the finitely long chains 
considered in the simulations have a helical ground state. 
The exact calculation of the ground state energy for 
chains of finite length is a formidable task. We have 
estimated an upper bound with the method of simulated 
annealing. At chain length TV = 50, we get U toi < — 210e, 
corresponding to U tot /N < — 4.2e per monomer, which 
is well below the energy of an optimized helix. The low 
energy conformations do not even exhibit local helical or- 
der. Hence we conclude that chains of length N = 50 or 
longer are not helical: The helical state observed in the 
simulations is a kinetic trap. 

Intriguingly, it does not even correspond to the en- 
ergetically optimized helix. The energy minimization 
for helical structures (minimized with respect to the ra- 
dius, the pitch, and the bond length) yields a series of 
local minima characterized by different numbers n of 
monomers per turn, i.e., n = 3.92, 4.81, 5.75, 6.71, •• •. 
The lowest minimum is at n — 4.81. The structures of 
the helices shown in Fig. Q] are closer to the structure 
corresponding to n = 3.92 with an energy per monomer 
U tot /N = -3.04e. 

Insight into the initial process of helix formation can 
be gained from a linear stability analysis of the starting 
conformation, the stretched chain. The unperturbed ref- 
erence state is an energetically relaxed straight chain with 



bond length b ~ bo. The Eigenmodes of the Hessian ma- 
trix are the longitudinal (L) and transverse (T) phonons 
of this chain. Neglecting chain end effects, they have the 
wavelengths A m = bN/ m and the frequencies 
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(coB(27r6/Am)-l), (6) 



where d 2 V IjJ /dr 2 is negative and dV LJ /dr is positive. Un- 
stable modes are characterized by imaginary frequencies. 
If the spring constant a is sufficiently large, the straight 
chain conformation is stable with respect to all longi- 
tudinal phonons. However, it is unstable with respect 
to transverse phonons with short wavelengths, A < 4b. 
Hence helical modes with less than n = 4 beads per 
turn are unstable. The most unstable phonon mode cor- 
responds to a zigzag mode. Indeed, weak zigzags are 
observed at the onset of the collapse process, but helix 
formation soon takes over, due to the fact that the he- 
lical state corresponds to a local energy minimum. We 
conclude that the helix formation is initially a driven pro- 
cess. This explains why the system picks the helical state 
with n = 3.92, rather than the optimized helical state at 
n = 4.81. 

The subsequent collapse is analyzed most conveniently 
in a slightly modified system where one chain end is fixed, 
corresponding, e.g., to a grafted chain. The helix forma- 
tion then proceeds from only one end (the free end) , and 
the defect in the middle of the chain is avoided. Follow- 



ing Kemp and Chen [15j], we define the two helical order 
parameters 
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where the Ui are unit vectors proportional to Ui oc 
(rj — rj_i) x (r.; + i — i-;), and N is the chain length. The 
parameter Ha characterizes the global helical order, and 
H2 the local order along the chain. 

Fig. [2K) shows a typical time evolution of the or- 
der parameter Ha for such a chain at the temperature 
T = 0.04e/fcs. One can distinguish two stages: In the 
beginning, Ha rises linearly up to a saturation value. This 
corresponds to a driven regime, where first the instability 
described above triggers the formation of a helix at the 
free chain end, and then this helix propagates into the 
chain. In the second stage, the chain stays perfectly he- 
lical for a while, until suddenly part of the helical order 
gets lost. To analyze this second stage, we have recorded 
the lifetimes of the helical state in 200 independent runs 
(different random numbers), starting from an identical 
helical chain conformation. The lifetime was defined as 
the time when the energy dropped below a given thresh- 
old E t , which was chosen below the initial energy Eq 
such that Eq — E t was of the order le, much larger than 
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FIG. 2: Time evolution of helical order for endgrafted chains 
with length N = 50 at temperatures around fcsT/e ~ 0.04. 
(a) Evolution of order parameter Hi for one specific run, with 
corresponding chain conformations (grafting point is repre- 
sented by a dark sphere), (b) Number N(t) of remaining 
purely helical conformations after the time t in a sample of 200 
initially helical conformations, at temperature fcsT/e = 0.04 
for two chain lengths as indicated. Dashed lines show fits to 
an exponential, N(t) ~ exp(— t/to) with escape time to. (c,d) 
Escape rates l/to in units of 1/r as a function of chain length 
for the temperature fcsT/e = 0.04 (c), and as a function of the 
inverse temperature for the chain length N = 50 (d). Dashed 
line in (c) shows the result of a fit to l/to oc exp( — AS/feaT), 
giving AE = 0.041 ± 0.002e. 



the typical energy fluctuations of an intact helix (roughly 
0.002 e per monomer). To ensure that the results did not 
depend on the threshold, we carried out the analysis for 
at least two different threshold values that differed by 
le, and checked that the results did not change. This 
was done for a range of chain lengths, N — 25 — 100, 
and temperatures, ksT/e = 0.04 — 0.045. In all cases, 
the distribution of the lifetimes was nicely exponential. 
Two examples are shown in Fig. Eb)- We conclude that 
the decay from the "perfect" helical state is a stochas- 
tic, rate-driven process. According to Kramer's rate the- 
ory [lfj], the height of the energy barrier, AE, can be 
estimated from the temperature dependence of the es- 
cape rate, l/to oc exp(AE /ksT). The fit of our data to 
this law gives AE = 0.041e (Fig. Ht). Hence perfect he- 
lices should persist for a long time at temperatures well 
below AE/kB, and decay rapidly at temperatures above 
AE/ks- This is indeed observed in the simulations. 

More insight into the nature of the escape process can 
be gained from looking at the dependence of the escape 
rate l/to on the chain length N. We assume that the 
escape is initiated by the nucleation of some defect. If this 
defect is localized at the end of the chain, l/to should not 
depend on N; if it is localized somewhere in the middle, 
l/to should be proportional to N. Fig. [Hi) shows that 




FIG. 3: Time evolution of the helical order parameter Hi (a) 
and the total energy Utot (b) in units of e at temperatures 
between fcsT/e = and fcsT/e = 0.3 for initially stretched 
chains, averaged over 180 independent runs. The distance 
between levels in the contour plots is 0.05 in (a) (lowest level 
0.05) and 20e in (b) (lowest level -180e). 

the escape rate increases with N , thus the escape defect 
forms in the middle of the chain rather than at the end. 
However, the dependence is not strictly linear; the true 
situation is more complicated than our simple argument 
suggests. 

At temperatures above T = 0.04e/fcs, the "perfect" 
helices decay rapidly, but the system still retains heli- 
cal order. To determine the temperature range where 
long-lived helices appear, be they perfect or imperfect, 
we have calculated the average time evolution of H4 for 
temperatures up to T = O/Se/ks. The result is shown 
in Fig. [3^). Helices are observed at temperatures up to 
T ~ 0.1 — 0.15e/ks- It is worth noting that in this heli- 
cal regime, the average total energy U to t is higher than at 
higher temperatures ( Fig. [3b)) - This proves once more 
that the helical state does not correspond to a true free 
energy minimum fll\ . 

Up to now, hydrodynamic effects mediated by the sol- 
vent were disregarded. To assess the influence of the 
latter, we have also carried out simulations of a dynami- 
cal model that includes hydrodynamic interactions via an 
appropriate mobility tensor Dij . Eq. |(4]) is then replaced 
by 

h = J2 D vh + Yl&rf3, (8) 

3 j 

where fj is distributed as before (Eq. ([3|) and £ fulfills 
= D. When determining the mobility tensor for 
our chains, we must account for the fact that they are 
endgrafted: The first bead is subject to a constraint 
force that ensures v\ = in Eq. |[8|). In the absence 
of noise, the total force on bead 1 is thus given by 
fi = -Y,&i[ D uT lD i7fh where Dg- is the mobil- 
ity tensor for free chains. Hence the effective mobility 
matrix between the remaining beads is 

Di^Dff-DfrlDtr^Dtr- (9) 

The free mobility tensor D f "' was approximated by the 
Rotne-Prager tensor [3, EH] with solvent viscosity f. 
The 'square root' tensor £ was determined by means of a 
Cholesky decomposition of D. 
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FIG. 4: Left: Time evolution of the local helical order param- 
eter H 2 (main plot) and the global order parameter Hi (inset) 
in an initially stretched chain of length N = 50 at tempera- 
ture UbT = O.Ole without (dashed lines) and with (straight 
lines) hydrodynamic interactions, averaged over 100 samples. 
Right: Examples of corresponding conformations. 



interlinked mineral pores [20j. In such environments, hy- 
drodynamic interactions are largely screened, and chains 
collapsing from stretched conformations may form nicely 
ordered helices. In an active system where polymers con- 
stantly attach and detach to other objects and are per- 
manently stretched and released, helices might thus have 
appeared perpetually for purely kinetic reasons. They 
could then have been stabilized a posteriori by suitable 
chemical modifications. 

S.A.S. thanks the MHESR, Algeria, and the DFG, Ger- 
many (SFB 613) for financial support during an extended 
visit to Bielefeld. The polymer conformations have been 
visualized using the open-source package VMD |2ll |. 



Fig. [4] shows the average time evolution of the helical 
order in such a chain at T = O.Ole/fcs- At this tempera- 
ture, simple Langevin motion yields highly ordered long- 
lived helices (cf. Fig. [J). In the presence of hydrodynamic 
interactions, the overall motion of the chain is much more 
cooperative, and the local propagation mecanism for he- 
lices does not work. As a result, the perfect global helical 
order disappears, only local order remains. We note in 
passing that the most unstable chain mode, the zigzag 
mode, rises to much higher amplitudes and persists much 
longer in the presence of hydrodynamics than without. 

A more detailed account of the hydrodynamic simu- 
lations shall be given elsewhere. Here, we just conclude 
that hydrodynamic interactions largely reduce the effect 
reported above. The "helical trap" is avoided in the pres- 
ence of hydrodynamics. This is consistent with previous 
simulations of chain collapse [1, 0], which also showed 
that hydrodynamic interactions tend to inhibit trapping 
during chain collapse. 

To summarize, we have observed the spontaneous for- 
mation of well-ordered, long-lived helices in Brownian 
dynamics simulations of simple flexible model polymers 
at low temperatures. The effect relies on three con- 
ditions: (i) The quality of the solvent must be suffi- 
ciently bad, i.e., the effective attractive interaction be- 
tween monomers must be well above the thermal energy 
ksT . (ii) The chain must be straight initially. The ini- 
tial conformation then has unstable helical modes, which 
trigger the helix formation, (iii) Hydrodynamic inter- 
actions should be screened. This is the case, e.g., in a 
porous or molecularly crowded environment. Otherwise, 
the local ordering process is disturbed by global cooper- 
ative processes. Remnant helicity (local helical order) is 
still observed in the presence of hydrodynamic interac- 
tions, but the global order is much reduced. 

We speculate that this effect might have contributed to 
establish helices as one of the most important structural 
elements in biopolymers. Recently, it has been hypothe- 
sized that the early molecular evolution of life may have 
taken place in strongly confined environments, e.g., thin 
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